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Abstract 

In this paper, I will discuss the Runge-Kutta method of solving simple 
linear and linearized non-linear differential equations. I start by stating why 
the Runge-Kutta method is ideal for solving simple linear differential equa- 
tions numerically in comparison to more elementary methods. I will then 
proceed to explain what steps the method actually carries out in solving the 
differential equation along with the matlab code I used to write a simple 
Runge-Kutta solver and the output of the code, given some basic differential 
equations. I will then end by discussing the pitfalls of this way of solving 
differential equations. 

1 Introduction 

Any nth order linear ordinary differential equation can be expressed in the 
form of a system of first order linear ODEs. Numerical approximation can be 
applied to effectively analyse differential equations that may not be easy to 
solve by traditional methods. Once any ODE has been expressed as a linear 
system, any of the methods of linear algebra can be applied on the system. 
The Runge-Kutta method is a far better method to use than the Euler or 
Improved Euler method in terms of computational resources and accuracy. 
All these methods use a fixed step size, but there are other methods that 
use a variable step size (though not neccessarily better in all circumstances). 
While essentially the Euler methods are simple Runge-Kutta methods, I, like 
has now become common, refer to the fourth-order Runge-Kutta method as 
the Runge-Kutta method. 

2 The Euler and Improved Euler methods 

For an initial value problem 

dy/dx = f(x,y),y(a) = y 

that has a unique solution y(x) on the closed interval [a, b] and given that 
y(x) has a continous second derivative on that interval, then there exists a 
constant C such that for the approximations yi to the actual values y(xi) 
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computed using the Euler's method with step size h > 0, 

\y n -y(xn)\ < Ch 
For n = 1,2,3, . . . ,k 
for the Euler method, and 

\y n -y(x n )\ < Ch 2 

for the improved Euler method. Here, h is the step size. This implies that 
the error is of order h in the Euler method and order h 2 in the improved Euler 
method. The proof can be found in the book, Ordinary Differential Equa- 
tions by G. Birkhoff and G.C. Rota. On the other hand, the Runge-Kutta 
method is a fourth-order method (Runge-Kutta methods can be modified 
into methods of other orders though). The Euler methods suffer from big 
local and cumulative errors. The improved Euler method and the Runge- 
Kutta method are predictor-corrector methods and are more accurate than 
the simple Euler method. 

3 The Runge-Kutta Method 

This method uses the simple fact that, for a given actual change in the out- 
put y, we can use the fundamental theorem of calculus to express the change 
in the form of an integral 

y(x n+1 ) - y(x n ) = y'(x) alx = f^ +h y'(x) dx 

We can then use Simpson's rule for numerical integration 

y(x n +i) - y(x n ) ~ | (y'(x n ) + ty'(x n + |) + y'(x n+1 )) 

y n +i « y n + | (y'(x n ) + 2y'(x n + §) + 2y'(x n + |) + y'(x n+1 )) 

I used matlab to write out some code that implements an algorithm that uses 
the Simpson's rule to approximate solutions for a given number of steps. Us- 
ing the code, one can solve any nth order linear ODE as long as you specify 
the initial conditions. I have shown that by attempting to approximate the 
solutions to a first order ODE and a second order ODE. The code is shown 
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below: 



function [T,X] = hrungekutta(t ,x,tl ,n,f ) 
T = t; 
X = x' ; 

h = (tl - t)/n; 
for i = l:n; 

kl = f(t,x); 

k2 = f(t + h/2, x + kl*h/2); 

k3 = f(t+h/2, x + k2*h/2); 

k4 = f(t+h, x + k3*h); 

x = x + h*(kl + 2*k2 + 2*k3 + k4)/6; 

t = t+h; 

T = [T;t] ; 

X = [X;x'] ; 

end 

This function takes the form f(t in uiai, [x ; yo\,tf ina i, steps, @ode) 

I implemented ODEs as functions. So it can take in any ODE as long as you 
express it in the form of a linear system and you specify its initial values. 
Note that in Matlab, to input a function to another function, you must use 
the @ sign in front of the input function's name. 

For example, the function ddxeqnegx is the ode x" = —x . 
In this case, you will input @ddxeqnegx . 

Example: [T,X] = hrungekutta(0,[0;l], 5, 10, ©ddxeqnegx) 

This solves x" = —x whose general solutions are of the form [sin(t), cos(t)] 
for t running from to 5 where the initial values are x — 0, y — 1, since you 
change x" = —x to the system of first order ODEs x' — y,y' — —x 

The code for the input ODE is: 

function X = ddxeqnegx (t ,x) 
X=x; 

X(l) = x(2); 
X(2) = -x(l); 
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A graph of the solution is shown in the next page: 

Cosine and Sine waves produced by solving x"=-x in 1 000 steps 
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An example of a first order ODE is [x,y]=hrungekutta(0,[0],l,10,@pidiff) 
which approximates the solution to ^ 



(l+x 2 ) 



with the initial value yo — 0. 
Note that the general solution to this differential equation is y = — , 4 -, , . 

° ~\. a arctanix) 

Therefore, y\ ~ 7r . Outputs for a ten step Runge-Kutta method are shown 
below. 

x = 





0.1000 
0.2000 
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0.3000 
0.4000 
0.5000 
0.6000 
0.7000 
0.8000 
0.9000 
1.0000 





0.3987 
0.7896 
1.1658 
1.5220 
1.8546 
2. 1617 
2 . 4429 
2.6990 
2.9313 
3. 1416 

With 1000 steps, I got n = 3.141592653589791 



4 Pitfalls in the Runge-Kutta method and 
other numerical methods 

There are a number of problems faced by the Runge-Kutta method. While I 
will not go into the details here, I will use an example equation to illustrate 
a problem that one could face. 



Consider the equation 

g = 5y-6e-*,y(0) = l 
Its general solution is 
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y {x) = e~ x + Ce 5x 



With the initial condition y(0) = 1, then C = 

A small error in the evaluation of the output leads to a non-zero value of 
C. Since the exponent is big, a very small error in approximation (which is 
almost inevitable) leads to a huge error in the output. 
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